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Abstract. Core-collapse supernovae produce an intense burst of electron antineutrinos in 
the few-tens-of-MeV range. Several Large Liquid Scintillator-based Detectors (LLSD) are cur¬ 
rently operated worldwide, being very effective for low energy antineutrino detection through 
the Inverse Beta Decay (IBD) process. In this article, we develop a procedure for the prompt 
extraction of the supernova location by revisiting the details of IBD kinematics over the broad 
energy range of supernova neutrinos. Combining all current scintillator-based detector, we 
show that one can locate a canonical supernova at 10 kpc with an accuracy of 45 degrees 
(68% C.L.). After the addition of the next generation of scintillator-based detectors, the 
accuracy could reach 12 degrees (68% C.L.), therefore reaching the performances of the large 
water Cerenkov neutrino detectors. We also discuss a possible improvement of the SuperNova 
Early Warning System (SNEWS) inter-experiment network with the implementation of a di¬ 
rectionality information in each experiment. Finally, we discuss the possibility to constrain 
the neutrino energy spectrum as well as the mass of the newly born neutron star with the 
LLSD data. 
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1 Introduction 

Neutrinos 1 have already been detected from the type II supernova SN1987A that appeared in 
the Large Magellanic Cloud (LMC) in 1987. At that time three detectors, Kamiokande II [1], 
IMB [2] and Baksan [3] detected twenty-nine neutrino events, hence opening the field of 
extrasolar neutrino astronomy (see Ref. [4] for a recent review). 

Since they leave the dense medium of the core-collapse supernova at the early stages of 
the process, neutrinos precede the optical photons by several hours [5]. An early pointing 
could allow the astronomical community to quickly locate the supernova and consequently 
study the explosion early phases with the world’s telescope armada. 

1 For simplicity reasons and unless stated otherwise, the term neutrino will be used for both electron 
neutrinos and electron antineutrinos in the following. 
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The estimated core-collapse supernova rate in our Galaxy is a few per century. The 
possible progenitor could be located anywhere, nevertheless the supernova probability dis¬ 
tribution computed from the red giant stars distribution is centered around 11.9 kpc with a 
RMS dispersion of 6 kpc [6]. 

Neutrino scattering on electrons in a very large water Cerenkov detector, such as Super- 
Kamiokande [7], provides the best way to locate a supernova through neutrinos, with an 
accuracy of about 8 degrees at 10 kpc (95% C.L.) [8]. However, the pointing accuracy is 
degraded by the IBD reactions, with a neutrino scattering signal over IBD noise ratio of 
about 1/30. On the other hand, the large cross section of the IBD reaction compared to 
neutrino scattering allows the detection of supernova v e in smaller detectors. 

In this work, we revisit the early pointing of a core-collapse supernova using the IBD 
reaction in LLSD. Although it does not provide an event by event pointing ability, the high 
cross section of this reaction leads to a high number of interactions, almost free of backgrounds. 
The statistical analysis of the angular distributions of the reaction products (positron and 
neutron) can be correlated to the neutrino direction and its associated uncertainty, enabling 
an early pointing of the supernova. 

While a single experiment can independently provide a most probable supernova loca¬ 
tion, we show that the prompt combination of data from all large liquid scintillator detectors 
could reduce the held of view for further optical-based astronomical observations. We con¬ 
clude by discussing other astrophysical information that can be retrieved such as the mass of 
the newly born neutron star and the supernova distance. 

2 Core collapse supernova neutrino signal 

In this section we outline the physics of core collapse supernovae and the associated neutrino 
signals. The present description is by no means complete and further details can be found in 
Ref. [9] for instance. 

Core-collapse supernovae are massive stars that undergo a gravitational collapse at the 
end of their life, emitting 99% of the gravitational energy through neutrinos and antineutrinos 
of all flavors. A supernova neutrino luminosity curve closely follows the different phases of 
the explosion, from the tens-of-milliseconds burst of v e (due to electron capture on nuclei), 
to the hundreds of milliseconds accretion phase and the ten-second-long cooling phase of 
the left proto-neutron star. During accretion and cooling, neutrinos and antineutrinos of all 
flavors are emitted with (E Ue ) < {E Ue ) < (E Ux ) with v x = i/ r , v T . This hierarchy of 
average energies reflects neutrino interactions with matter at tree level, and the fact that the 
star is neutron-rich. At the neutrinosphere boundary, where neutrinos start free-streaming, 
the neutrino fluxes are well described by a Fermi-Dirac or a power-law distribution [10]. 
Outside the neutrinosphere, the Mikheev-Smirnov-Wolfenstein (MSW) effect [11, 12] is well 
established. Additional flavor conversion phenomena occur due to neutrino self-interaction 
as well as the presence of shock waves and turbulence. While many features concerning the 
mechanisms and the occurrence of such flavor conversion phenomena have been unraveled, 
important open questions remain [13]. Note also that the recently discovered lepton-number 
emission self-sustained asymmetry (LESA) is expected to produce a dipole pattern during 
the accretion phase [14], In fact, as we will show, the shape of the neutrino spectra has a 
negligible impact on the directional reconstruction of a core-collapse supernova. 

In what follows we focus on the v e emission, relevant for IBD detection. We are consider¬ 
ing three examples of neutrino flux prediction: a maxwellian-like distribution, the Livermore 
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flux [15], and the GVKM (Gava-Kneller-Volpe-McLaughlin) flux [16]. Figure 1 displays the 
three expected energy spectra. We will show that the shape of these distributions only has a 
negligible effect on our following results. 

The maxwellian distribution is a common thermal distribution used to simply describe 
the energy spectrum of the neutrinos immediately after the burst. The Livermore distribution, 
similar to a Fermi-Dirac or a power-law distribution, describes the v e energy spectrum after 
exiting the supernova neutrinosphere, without conversion or oscillation effects. The GVKM 
distribution provides a modeling of flavor conversion due to neutrino self-interactions and 
shockwave effects that affect the fluence and the shape of the v e spectrum. Details of these 
effects can be found in [16]. It reduces the number of expected electronic antineutrino events 
by 50% with respect to the Livermore flux. It is an improved calculation of the expected 
energy spectrum and we therefore consider it as our reference flux. 


xIO 9 



Figure 1. Expected core-collapse supernova D e emission spectra (integrated over 10 s) based on the 
Maxwell, Livermore, and GKVM fluxes at 10 kpc. The three spectra are normalized to the same 
energy released by the supernova, 3 x 10 52 ergs (2 x 10 58 MeV), equally distributed among each 
neutrino and antineutrino species. The integrated fluences of the three distributions are respectively 
9.7 x 10 10 i> e .cm -2 , 1.9 x 10 11 P e .cm -2 and 9.7 x 10 10 v e .cm~ 2 . 


3 Neutrino detection in large liquid scintillator detectors 
3.1 Neutrino detection through IBD 

Core-collapse supernova neutrinos can interact through charged-current reactions on nuclear 
targets. Here we focus on the IBD in which only v e participate, since it leads to the highest 
number of events in liquid scintillator detectors. Electronic antineutrinos are detected via 
IBD on free protons: v e + p —> e + + n for Ep e > E t ^ r ~ 1.806 MeV in the laboratory frame 
where the free proton is at rest. 

The IBD reaction cross section can be precisely computed from the V-A theory of weak 
interaction. However, depending on the considered energy scale, some approximations are 
often performed to simplify the mathematical expressions. For core-collapse supernova v e 
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energy range, up to Ep e —100 MeV, we consider the cross section computed by A. Strumia 
and F. Vissani in Ref. [17] and displayed in Figure 2: 


da G\ cos 2 9 C . |2 

d t 2TT(s-m 2 p ) 2W 1 ’ 


(3.1) 


where Gf is the Fermi constant, cos 9c is the Cabbibo angle, m p is the proton mass, A4 is 
the invariant amplitude for the IBD reaction explicitly given in Ref. [17], s = (p u + p p ) 2 and 
t = (jpv — p e ) 2 are expressed as functions of the neutrino, proton and positron 4-momenta p u , 
p p and p e . Other cross section expressions such as a naive “p e E e ” relation (see Ref. [17]) and 
computations carried out by P. Vogel and J. Beacom [18] are displayed in Figure 2 as well. 
As stated in Ref. [18], the total cross section computed by Vogel-Beacom is only valid up to 
60 MeV. 



Figure 2. Left: IBD cross section based on three different approximations that differ significantly for 
energies greater than 30 MeV. Right: Expected neutrino spectrum for a supernova located at 10 kpc, 
normalized to the interaction on one free proton. 


Following the computation of Ref. [17] the outgoing positron and incoming v e energies 
are related by: 

( E u — h) (1 + e) + e cos 9 e \J(E u — h) 2 — m^K 


E e = 


K 


(3.2) 


with 


5 = 

(m\ — m p — ml) /2 m p 

(3.3) 

e = 

Eulm p 

(3.4) 

K = 

(1 + e) 2 - (ecos 8 e ) 2 , 

(3.5) 


where E v is the neutrino energy, E e is the outcoming positron total energy, 0 e is the angle be¬ 
tween the positron and the incoming antineutrino directions in the laboratory frame and m n , 
m p and m e are respectively the neutron, proton and electron masses. We call visible energy 
the positron energy deposited in the detector by ionization and annihilation, corresponding 
to -Fyis — E e -{- 77l e . 
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After its emission, the neutron undergoes several elastic scatterings in the hydrogen-rich 
medium of the liquid scintillator. This process called moderation decreases its kinetic energy 
until it reaches typical thermal equilibrium at ~0.025 eV. After thermahzation, the neutron 
diffuses before its capture on a nucleus. The average duration of this diffusion process depends 
on the liquid composition. It varies from tens of microseconds for a Gadolinium(Gd)-doped 
liquid to hundreds of microseconds for a non-doped liquid. The use of Gd doping increases the 
detection efficiency while reducing the impact of accidental background. However, it is worth 
noting that despite the Gd-doping, a fraction of the neutrons will be captured on hydrogen 
(H), depending on the liquid composition. The current technologies developed for the recent 
reactor neutrino experiments [19-21] guarantee the possibility of doping the next generation 
of large detectors with Gd at the multikiloton-scale. 

An IBD event in a LLSD is thus characterized by a prompt positron event which deposits 
a visible energy between 1 and up to 100 MeV, followed by a delayed burst of 7 -rays arising 
from neutron capture on a target nucleus within r ~ 10 — 1000 fjs, depending of the target 
liquid composition. These prompt and delayed events are spatially correlated, within a few 
cubic meters. 



Figure 3. Left: IBD interaction rate in various LLSD as a function of the visible energy, for a core 
collapse supernova located at 10 kpc. “Future only” and “All” are respectively the combination of 
LENA and JUNO and the combination of all LLSD, current and future. Right: IBD interaction 
rate for LLSD as a function of the supernova distance. The number of expected IBD events in 
Super-Kamiokande is shown for comparison. 


3.2 Running detectors 

The detection of supernova neutrinos in liquid scintillators relies on the existence of numerous 
large scale detectors with a mass bigger than several tens of tons. These volumes are large 
enough to detect a significant number of events and tag a supernova in our Galaxy or its closest 
satellites, such as the small and large Magellanic Clouds, with a negligible fake detection rate, 
thanks to the coincidence between detectors. The non-exhaustive list of detectors considered 
in this work and their relevant characteristics are collected in Table 1. Their respective 
expected event rates are displayed on Figure 3 and Table 2. 

NOi^A (NuMI Off-Axis v e Appearance) is an experiment dedicated to observe v e ap¬ 
pearance using Fermilab’s 17 beam NuMI [22]. It operates above ground with two detectors, 
each made of individual cells of liquid scintillator of 3.87 cm x 6.00 cm x 15.7 m. With 
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Table 1 . Main detector features used in our simulation. 


Detector 

Dimensions (m) 

Shape 

V (nP) 

free H/nr 5 

# free H 

LVD 

a=10, b=13.2, c=22.7 

Box 

1020 

7.46 ,10 28 

7.6 ,10 3i 

KamLAND 

r=6.50 

Sphere 

1150 

6.60 .10 28 

7.6 .10 31 

SNO+ 

r=6.00 

Sphere 

904 

6.24 ,10 28 

5.6 ,10 31 

Borexino 

r=4.25 

Sphere 

321 

5.30 .10 28 

1.7 .10 31 

Daya Bay 

r=2.00, h=4.00 

Cylinder 

50(x6) 

6.24 ,10 28 

1.8 .10 31 

Double Chooz 

r=1.70, h=3.55 

Cylinder 

32(x2) 

6.55 .10 28 

0.4 ,10 31 

Reno 

r=2.00, li=4.40 

Cylinder 

55(x2) 

6.24 ,10 28 

0.7 .10 31 

NOvA (far) 

a=15.7, b=15.7,c=132 

Box 

32500 

6.24 ,10 28 

2.0 .10 33 

NOvA (near) 

a=3.5, b=4.8, c=9.58 

Box 

161 

6.24 ,10 28 

1.0 .10 31 

MiniBooNE 

r=5 

Sphere 

523 

7.31 .10 28 

3.8 .10 31 


12,480 cells, NOi'A’s near detector consists of 145 tons of active material while the far detec¬ 
tor and its 761’856 cells encloses 24 ktons of liquid scintillator. Its important background rate 
at sub-GeV energies and its energy threshold, above the neutron capture on H at 2.2 MeV, 
does not allow an efficient neutrino burst detection nor a neutron vertex reconstruction yet. 
NOizA will thus not be taken into account in the following direction reconstruction analyses. 

The Large Volume Detector (LVD)[23] is located at the Laboratori Nazionali del Gran 
Sasso (42° 27’ 10" N, 13° 34’ 30" E), at the depth of 3500 m.w.e.. It is a 1 kt liquid scintillator 
detector in the form of an array of 840 scintillator counters, 1.5 m 3 each. While its active 
mass allows the detection of an important number of neutrinos, the size of LVD’s individual 
scintillating counters prevents any precise position reconstruction. Hence, like NOz'A, it will 
not be taken into account in the following directionality analyses. 

From now on, the following detectors will be considered in our directionality analyses. 
The KamLAND experiment [24] is a 1 kton liquid scintillator (80% dodecane and 20% pseu¬ 
documene) balloon located within the Kamioka mine in Japan (36° 25’ N, 137° 18’ E), at 
2700 mwe under the Ikenoyama mountain. 

The SNO+ experiment [25] is currently in its final construction phase at the SNOLAB 
laboratory (46° 28’ 30" N, 81° 12’ 04"), shielded with an overburden of 6000 m.w.e.. SNO+ 
will be filled with 780 tons of Linear Alkyl Benzene (LAB). 

The Borexino experiment [26] is located at the Laboratori Nazionali del Gran Sasso in 
Italy (42° 27’ 10" N, 13° 34’ 30" E) and shielded with an overburden of 3500 m.w.e.. It 
consists of 278 tons of pseudocumene-based liquid scintillator. 

Three reactor neutrino oscillation experiments, Daya Bay [19], Double Chooz [20], and 
RENO [21] are currently taking data to measure the #13 mixing angle. They all consists 
of several identical detectors located within two kilometers from several powerful nuclear 
power plant cores used as intense electron antineutrino sources. The Daya Bay experiment is 
located in the Guang-Dong Province, on the site of the Daya Bay nuclear power station (22° 
36’ 59" N, 114° 32’ 28.1" E). Each of the six cylindrical identical detector modules contains 
an effective volume of 20 tons of 1 g/1 Gd-loaded liquid scintillator, and 22 tons of non-doped 
scintillator. The Double Chooz experiment is located close to the twin reactor cores of the 
Chooz nuclear power station located in the French Ardennes (50° 08’ 43" N, 4° 80’ 165" E). 
The two cylindrical identical detectors contain a 8 -ton fiducial volume of liquid scintillator 
doped with 1 g/1 of Gd and 18 tons of non-doped scintillator. The RENO experiment is 
located on the site of the Yonggwang nuclear power plant in Korea, about 400 km south of 
Seoul (35° 24’ 04" N, 126° 25’ 31.5" E). The two cylindrical identical detectors contain a 
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16-ton fiducial volume of liquid scintillator doped with 1 g/1 of Gd and 39 tons of non-doped 
scintillator. Given the large volume fraction of undoped scintillator liquid in those three 
detectors, in the following we only consider the events expected from neutron capture on H. 
Considering the energy and vertex reconstruction enhancement brought by neutron capture 
on Gd, this is a conservative approach. 

The MiniBooNE detector [27], located at Fermilab (41° 83’ N, - 88 ° 26’ E), is a spherical 
volume containing 680 tons of mineral oil built at ground level. The project of adding PPO 
in the mineral oil proposed in Ref. [28] indicates the detector might still be running for a few 
years. The addition of PPO would enhance scintillation and set MiniBooNE in our LLSD 
list. 

Table 2. Expected number of interactions taking into account the exact liquid scintillator chemical 
compositions and densities, for a 10 kpc core-collapse supernova. 



Livermore 

GVKM 

LVD 

319 

190 

KamLAND 

318 

189 

SNO+ 

237 

141 

Borexino 

71 

43 

Daya Bay 

79 

47 

Double Chooz 

18 

11 

RENO 

29 

17 

NO^A far 

6811 

4057 

NO^A near 

32 

19 

MiniBooNE 

161 

96 


3.3 Forthcoming detectors 

Beyond the currently running detectors two very large liquid scintillator detectors projects 
are being considered: JUNO and LENA. These detectors, besides providing an answer to 
some of the major open questions remaining in the neutrino held, are perfect apparatus for 
studying core-collapse supernovae. 

JUNO (Jiangmen Underground Neutrino Observatory) is under construction and will 
consist of a 20 kt LAB-based liquid scintillator detector located at Jiangmen (22°37’ N, 
112°70’ E) in South China [29]. Being 20 times as massive as SNO+, the same factor on the 
expected event number is expected. 

LENA (Low Energy Neutrino Astronomy) is a project designed to be a next-generation 
liquid scintillator detector on the scale of 50 kt [30]. Its location is still under investigation 
but the most probable is the Pyhasalmi mine (63°39’ N, 26°2’ E, 4000 m.w.e.), located in the 
Pyhajarvi province in Finland. 

4 Supernova early pointing in LLSD 

Our goal is to determine the supernova direction in the sky from a statistical analysis of 
all interactions associated to the electron antineutrino burst, taking advantage of the IBD 
kinematic features keeping some memory of the incoming neutrino direction. For v e IBD 
interactions in LLSDs, prompt and delayed event features, such as energy and interaction 
vertices, are determined on an event-by-event basis in each detector. For each event, one can 
then define a vector starting at the prompt vertex and ending at the delayed vertex. This 
information is complemented by the prompt visible energy. 
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4.1 Relevant IBD kinematics features 

In this section we review the IBD relevant kinematic features used to constraint the v e in¬ 
coming direction. We compare analytical calculations of Ref. [18] with our toy Monte Carlo 
simulation described in Appendix A. 

Visible energy 

The relation between the positron energy and the neutrino energy is straightforward at reactor 
neutrino energies, E v <10 MeV, but it gets ambiguous at higher energies. Eq. 3.2 shows that 
E v both depends on E e and cos 9 e . Since only E e can be measured, a degeneracy appears 
in the v e energy estimate. For reactor neutrinos (E u <10 MeV), this degeneracy is weak, as 
shown by Figure 4, while it strengthens at higher energies. Since cos 9 e cannot be measured, 
this induces an error of up to 10% in the v e energy reconstruction. The positron takes 
the maximum of the energy when it is emitted along the neutrino direction. On the other 
hand, it goes preferably slightly backwards when taking a relatively small energy, the rest 
being available to the neutron that can therefore scatter further forwards with larger emission 
angles. 
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Figure 4. Relationship between E v and E e as a function of the positron emission angle, cos 6 e . The 
grey bands correspond to the range of possible reconstructed neutrino energies for a given positron 
energy. 


Positron and neutron directions 

The outgoing positron angular distribution is obtained by multiplying Eq. 3.1 with the fol¬ 
lowing Jacobian: 

df __ 2 p e E v _ 

d cos 0 e m p + E v (1 — E e /p e cos 9 e ) 

For each pair (E u ,cos9 e ), E e is obtained from Eq. 3.2. For E u < 15 MeV, the positron 
is emitted on average in a slightly backward direction, whereas it moves preferably forward 
as the incident neutrino energy is higher, approximately following the form < cos 9 e >~ 


(4.1) 
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—0.034 \J— ml jE e +2.4 E v /m v computed in Ref. [18]. The neutron is always emitted in 
the forward hemisphere, with a maximum allowed angle increasing with the neutrino energy. 
Figure 5 shows the mean values of cos0 e and cos 8 n obtained from our toy Monte Carlo 
simulation compared with the Vogel-Beacom analytic computation of the differential cross 
section, valid up to ~ 150 MeV. The positron-neutron vector, our directionality observable, 
is thus energy dependent. 




Figure 5. Left: Average cosine of the outgoing positron angle, cos 9 e , as a function of neutrino energies 
between E t h and 100 MeV. Right: Average cosine of the outgoing neutron angle, cos 9 n . Both panels 
show results of our toy Monte Carlo based on the Strumia-Vissani cross section compared to the 
analytical evaluations of Vogel-Beacom. 


Positron and neutron mean paths 

After their creation both particles move into the detector with respect to the interaction ver¬ 
tex before being annihilated or captured. 

As it crosses liquid scintillator, the positron loses its energy by ionization until it annihi¬ 
lates with an electron thus creating a pair of gamma rays. The mean path length between the 
positron creation at the interaction point and its annihilation obtained from our simulation 
is displayed in Figure 6, for different positron energies. As explained in Appendix A, the 
path length must be corrected to account for the detector reconstruction since the track of 
deposited energy is reconstructed as a point-like vertex. 

Naturally the mean path length increases with the energy, with a peak value ranging 
from a few tens of mm to a several hundreds of mm at higher energies. Figure 7 shows our 
results for the mean positron path length projected onto the incoming neutrino direction as 
a function of the neutrino energy. It ranges from a few mm for E v <15 MeV up to about 
40 mm at 80 MeV. Note that for low energy neutrinos (E u <10 MeV), the projected path 
length is ~ 0.05 mm, in agreement with the expected value for reactor neutrinos [31]. 

The neutron, emitted simultaneously with the positron, loses its energy through elastic 
scattering off light nuclei of the scintillator medium until it thermalizes with the medium. This 
non-isotropic process, depending both on the neutron energy and on the liquid composition, 
preserves a memory of the initial neutrino direction. The mean path length between the 
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Figure 6. Probability density function of the path length of positrons with different kinetic energies. 



Figure 7. Mean positron and neutron paths projections along the v e direction as a function of the 
neutrino energy. 


neutron creation and its thermalization obtained from our simulation is displayed in Figure 8 
(left), for different neutron kinetic energies. After being moderated the neutron scatters at 
thermal energies until it gets captured. The path length of this diffusion process can be seen 
in Figure 8 (right). 
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Figure 8. Left: Mean neutron path length before diffusion as a function of the neutron kinetic energy 
T n . Right: Diffusion length of thermalized neutrons in liquid scintillator with and without Gd-doping. 


The diffusion is an isotropic random walk spreading the distribution of the neutron cap¬ 
ture point. Our simulations have shown that the variation in path length from a non-doped 
scintillator liquid compared to another is negligible. Figure 7 shows our results for the aver¬ 
age neutron path length projected onto the incoming neutrino direction as a function of the 
neutrino energy. It ranges from 25 mm for E u <15 MeV up to about 80 mm at 80 MeV. It is 
worth noting that we recover the results commonly used for reactor neutrino experiment for 
which the projection is expected to be ~20 mm for neutrino energies up to 10 MeV. Most of 
the neutron captures occur on an hydrogen atom thus emitting a 2.2 MeV gamma ray whose 
energy deposition constitutes the delayed event. Since this gamma is emitted isotropically 
and its mean path length in the liquid is about 20 cm, it smears the capture vertex location. 


4.2 Direction retrieval 


Method 


The information returned by our simulation is only the visible energy and the reconstructed 
vertices of the prompt and delayed events, likewise a real supernova observation in a LLSD. 
From these vertices, we obtain the Delayed-Prompt unit vector X for each neutrino event 
having its origin at the positron reconstructed position and pointing to the neutron capture 
reconstructed location, as displayed in Figure 9. 

We then compute a direction ((f), 0) for each event with <f> and 6 respectively the azimuthal 
and zenith angles in the detector reference frame. The angular distributions extracted from 
a supernova simulation are displayed in Figure 10. These distributions are respectively the 
azimuthal angle distribution, the zenith angle distribution and the cosine projection on the 
neutrino axis expressed as: 


cos <5 


X ■ V 
VII x ||z7| 


(4.2) 


with cosd the normalized scalar product of X the Delayed-Prompt vector and V the neu¬ 
trino direction. All simulations were carried out using a conservative vertex resolution of 
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Figure 9. Sketch of the Delayed-Prompt vector X along the direction of the incoming neutrino. 


20 cm. Tests with a 15 cm vertex resolution were performed without yielding significant 
improvements. 


180 
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Figure 10. Probability density functions of azimuthal angle (top), zenith angle (center) and cosine 
projection on the neutrino direction (bottom) distributions extracted for a supernova pointing towards 
the (0,0) spherical coordinates in the detector frame. The angular distributions are centered on the 
(0,0) spherical coordinates. Neutrino energies are drawn from a GVKM spectrum. 

In order to consider all the relevant features of the IBD kinematics we reconstruct the 
direction with a x 2 minimization method comparing a model and a supernova dataset. Both 
are simulations of supernova neutrinos, performed with our toy Monte Carlo described in 
Appendix A, with a user-defined energy spectrum and cross section. The model dataset is 
computed with a large number of events in order to smoothen statistical effects and gener¬ 
ated along an arbitrary direction. On the other hand, the supernova dataset has a realistic 
statistic corresponding to the number of IBD interactions in a given detector and is generated 
along a user-defined direction (<j) q, 9q). By rotating the model distribution across the sky, our 
algorithm finds the best-fit angles (<f>, 6) corresponding the minimum value of the x 2 expressed 
as: 
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(4.3) 


N N 

X 2 (</>, e ) = £ £ (Yi - M. i (</>, e )) {v ~ 1 ) i . (Y) - m,- (0,0)) , 

* j 

where Y) and Mi (q i, 0) are the concatenation over N bins of the three angular distributions 
(azimuthal, zenith and cosine projection) of the supernova and model datasets, respectively. 
Vij is the N x N covariance matrix made of two independent blocks accounting for the bin-to- 
bin statistical uncertainty of the azimuthal and zenith angle distribution and a third block, 
correlated to the first two, accounting for the bin-to-bin statistical and systematic uncertain¬ 
ties of the cosine projection distribution expressed in Eq. 4.2. The V t j covariance matrix has 
been generated using 2,000 independent datasets of 5,000 neutrino events each, then normal¬ 
ized to one event. The resulting unitary covariance matrix is then normalized to the number 
of events in the supernova dataset being used in each y 2 minimization. 

The best-fit angles ((/), 6) corresponding to a minimum of x 2 are then stored and the 
process is repeated with a new supernova dataset having the same statistic. From all these 
trials, we obtain distributions of reconstructed azimuthal and zenith angles whose means and 
variances correspond to the reconstructed supernova directions and associated angular errors. 

Results with a single detector 

We now apply the method previously described to a single detector. We considered Kam- 
LAND in this first application. We generated several datasets with different number of events 
and computed their associated angular uncertainty using the method previously described. 
Figure 11 shows the interpolation of our results, best fitted by an inverse square root function, 
and describes the behavior of the angular uncertainty with respect to the number of detected 
events. Using this curve and the expected number of events computed from Section 3, we 
associate an angular uncertainty on the supernova reconstructed direction for KamLAND at 
different distances as shown in Figure 11. Several simulations were carried out using the 
different flux shapes presented in Section 2 and no sensible differences in the results have 
been found. With a number of events leading to an angular uncertainty greater than 90°, 
we consider that no directional information can be extracted and consequently discard the 
concerned detector from our directionality reconstruction. On the other hand, the detection 
of more than 50,000 events would lead to a trigger rate too high to be withstood, hence we 
arbitrarily consider it as our detection limit. Therefore, any observation of a higher neutrino 
rate in a single detector would, in our analyses, not bring additional direction information 
and the reconstruction error could not reach lower values. Between these two limits, the 
angular uncertainty increases with the distance and with the inverse of the number of events. 
The statistical nature of the IBD directionality process is well shown since several thousands 
of events are required to obtain a decent reconstruction uncertainty. As a consequence, the 
pointing potential of a single detector such as KamLAND gets degraded as soon as the disc- 
tance exceeds a few kpc, hence the interest of combining several detectors to improve the 
accuracy. 

Combining large liquid scintillator detector data 

Higher sensitivities on supernova position reconstruction are reachable adding the directional 
information provided by several detectors. Each of the detectors previously described in 
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Figure 11. a: Angular uncertainty as a function of the number of detected events, b: Angular 
uncertainty as a function of the distance in kpc in the KamLAND detector. 


Section 3 could detect a supernova assuming it leads to a sufficient amount of detected events. 
The addition of these signals increases the supernova detection potential in two ways. 

First, observing a coincident supernova signal in several detectors across the globe 
strengthens the confidence that a supernova has occurred. This is the current role of the 
SNEWS network that will be discussed in Section 5. Furthermore the time difference be¬ 
tween the supernova signals coming from those detectors could give us a first hint of localiza¬ 
tion using the triangulation method as discussed in [32], This method might not be accurate 
enough at typical supernova distances since it requires a good understanding of the luminosity 
time distribution and is strongly statistic dependent. A significant work of time calibrations 
between detectors would also be mandatory. 

Second, one can combine the reconstructed angular directions and errors provided by 
those detectors to improve the localization of the detected supernova. 


Table 3. List of individual detectors and their corresponding datasets. 



Existing 

SNEWS 

Near future 

All 

LVD 





KamLAND 

X 

X 

X 

X 

SNO+ 

X 

X 

X 

X 

Borexino 

X 

X 

X 

X 

Daya Bay 

X 

X 

X 

X 

Double Chooz 

X 


X 

X 

RENO 

X 


X 

X 

NOvA (far) 





NOvA (near) 





MiniBooNE 



X 

X 

JUNO 



X 

X 

LENA 




X 


In the context of a global analysis with N detectors combined, all the individual angular 
“cones” of lcr azimuthal and zenith uncertainties must be combined into a single one. To do 
so, we combined the individual angular errors and expressed the average angular uncertainty 
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as: 


a = 



(4.4) 


with ai the la uncertainty on the supernova direction reconstructed in the i th considered 
detector. 

In the following analysis, we considered 5 sets of combined detectors, summed up in 
Table 3. The “Existing” dataset is the combination of all current LLSD, “SNEWS” consists 
of the four LLSD belonging to the SNEWS network as detailed in Section 5, “Near future” is 
the sum of all current LLSD with the addition of MiniBooNE and JUNO, planned to operate 
within the next 5 years and finally “All” consists of the sum of all detectors, current and 
future. Results are displayed in Figure 12 for individual detectors and Figure 13 for the sets 
of combined detectors. 



-KamLAND (Ikt) 
-<3-SN0+ (0.8kt) 

- Borexino (0.3kt) 

—I— DayaBay (0.3kt) 

- □ - DoubleChooz (0.05kt) 
-A- RENO (0.1 kt) 

- 0- MiniBoone (0.7kt) 

- x - JUNO (20kt) 
-V-LENA (50kt) 


Figure 12. Angular uncertainty as a function of the supernova distance for different detectors and 
their associated masses. 


Considering a supernova located at 10 kpc, no current detector could provide a valuable 
direction information other than the hemisphere in which the supernova occurred. However, 
when combining several of these detectors, this same supernova direction is likely to be re¬ 
constructed within a 45°(68% C.L.) cone. While this angular opening is still larger than most 
of the optical instruments can accommodate, the addition of JUNO in a near future would 
reduce the uncertainty to 12°(68% C.L.), hence reaching a relevant level, almost competi¬ 
tive with the current Super-Kamiokande detector. With the conservative assumption that 
any galactic gravitational stellar collapse will occur within a 20 kpc radius, the addition of 
current and future detectors could reconstruct its position within a 14°(68% C.L.) cone. 
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SNEWS 
-© - Existing 
- - - Near future 
All 


Figure 13. Angular uncertainty as a function of the supernova distance for different detector com¬ 
binations. The mean reconstruction error for combined detectors is computed using the combination 
of the angular reconstruction error associated with each individual detector. Due to the small event 
rate expected in Double Chooz and RENO, the “SNEWS” and “Existing” datasets yield very similar 
results. 


5 Improving the SNEWS network? 

The coincident observation of a burst in several neutrino detectors would be a robust early 
warning of a forthcoming visible galactic core-collapse supernova. The main goal of the 
SuperNova Early Warning System (SNEWS) is to provide the astronomical community with 
this early alert of a supernova [33]. Each experiment belonging to SNEWS has its own 
supernova real-time trigger system dedicated to detect a burst of events that might be linked 
to a supernova. When the trigger requirements are met, an alert is sent to the global SNEWS 
server and is compared with possible alerts coming from other experiments. If a coincidence 
is found within a time window of 10 seconds, a global alert is sent to the SNEWS mailing 
list, universally accessible via snews.bnl.gov. In 2014, Daya Bay joined the SNEWS network 
already composed of Super-Kamiokande, LVD, Borexino, SNO+, KamLAND and IceCube. 
Adding Double Chooz, RENO and NOzzA to the network would increase the number of 
detected supernova neutrinos, thus strengthening the confidence level of an hypothetical alert. 

The SNEWS network is set up to provide information that a supernova has been detected 
via its neutrinos to astronomers. Although this information is mostly characterized by the 
observation of a burst of events in coincidence in several experiments, a directional information 
can be sent as well, if provided. Without this information, the supernova location remains 
unknown until its observation. Knowing this position, even roughly, over the sky would 
give astronomers a consequent head start and could lead to the observation of a complete 
supernova process, including luminosity rise. 

As shown in Section 4.2, even though current LLSDs in SNEWS can provide directional 
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information for nearby supernovae, Super-Kamiokande remains the best detector to locate 
them. However, with the replacement of SNO heavy water with liquid scintillator and before 
the construction of Hyper-Kamiokande [34], it stands alone as the world’s only large Cerenkov 
detector. Relying solely on Super-Kamiokande might not be enough to ensure a permanent 
watch on supernovae. Any maintenance operation or reconstruction software issue could put 
astronomers in the blind in case of a long-awaited supernova detection. 

With the help of SNEWS’s current and future liquid scintillator detectors, this scenario 
is unlikely to happen as the directional information of several detectors that observed the 
supernova can be combined. In the case of all SNEWS detectors operating at once, the 
Super-Kamiokande direction reconstruction will be improved by the combined reconstruction 
of all LLSDs. With only 4 LLSDs (KamLAND, Borexino, SNOj- and Daya Bay) capable 
of direction reconstruction, SNEWS is currently able to locate a 8 kpc supernova within a 
barely practical 35°(68% C.L.) wide cone. However, with the addition of all current LLSDs 
presented in Section 3.2 as well as JUNO and MiniBooNE in a near future, the opening of 
this cone could be reduced to 9°(68% C.L.) thus reducing by a factor 15 the area of the sky 
region of interest. 

In order to join the SNEWS network, each experiment would have to meet a few re¬ 
quirements, one of the most important being the robustness of its data acquisition system 
(DAQ), i.e. its capability to quickly and efficiently tag a supernova burst. In order to send 
a reliable alert to the SNEWS global server, each detector’s DAQ needs to be modified to 
operate on real-time mode and detect any unusual trigger rate increase. In the case of NOi'A, 
the detector’s trigger is directly synchronized with the NuMI beam trigger hence the neces¬ 
sity of operating the detectors in a “free” mode without any external trigger as presented in 
Ref. [22], An example of the implementation of such a trigger for the Daya Bay detectors is 
detailed in Ref. [35]. 

Implementing a combined pointing ability within the network might be more challenging. 
While the detection of a supernova burst solely relies on the observation of a burst of events 
in several detectors in coincidence, its pointing towards the sky requires energy and position 
reconstructions of these events. Most detectors operate on a run-by-run basis and usually 
analyze large datasets after long periods of data taking. Providing a neutrino direction within 
a time scale of a few hours would require offline and real time data analysis. This might be 
achieved by installing separate and independent online data acquisition and analysis chains 
capable of supporting large data rates and fast energy and position reconstruction based on 
the current analysis. 

6 Retrieving astrophysical information from supernova neutrinos 

Beyond our main goal of determination of the supernova direction, we discuss in this section 
additional information that could be retrieved by LLSDs, using the IDB process. In particular, 
we consider the determination of the supernova distance and then the assessment of the 
neutron star mass. The latter was investigated in Ref. [36] in the case of SN1987A. Clearly 
it would also be of interest to determine how well one could identify the average energy and 
pinching parameters of the supernova neutrino fluxes associated with the different phases of 
the explosion (neutronization burst, accretion and cooling). This has been studied e.g. in [37] 
by combining one- and two-neutron emission channels in lead-based detectors. However a 
similar study for scintillator detectors goes beyond the scope of the present manuscript. Just 
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as an example, we show how well one could discriminate the two supernova neutrino fluxes 
used in the present study. 

6.1 Measurement of the supernova distance 

Since, for a given supernova, the number of detected events varies with the inverse of the 
squared distance, the determination of its distance sounds straightforward in principle. To 
better understand the impact of statistical and systematic uncertainties, we developed a 
Monte-Carlo simulation to generate and fit detected energy spectra with a x 2 minimization 
method. This allows us to compare fake data fitted with the GVKM flux under various 
hypotheses. The main systematics considered are uncertainties on the neutron star features. 
Indeed the number of emitted neutrinos depends on the difference between the progenitor 
star and the neutron star gravitational binding energies. However, the first term can be 
considered negligible and the neutrino fluence therefore depends on the neutron star binding 
energy 2 Eb = f 7 ^ with M ns and R ns , being the mass and radius of the neutron star [38]. 
In the literature these parameters are usually set to M ns = 1.4 M 0 and R ns = 10 km [39]. 
From the average observed neutron star masses and radii computed in Ref. [40] and Ref. [41], 
we conservatively consider a relative uncorrelated systematic uncertainty of 10 % for both 
parameters. We simulated GVKM-based supernovae bursts located from 1 to 10 kpc away 
and reconstructed their distances and associated uncertainties (stat. + syst.). Hereafter, a 
conservative energy resolution of 6 % in the whole energy range is taken into account. The 
distance is reconstructed without bias and with an uncertainty of 40% and 15% at 10 kpc for 
a 1 kt and a 20 kt detector, respectively. These results are valid for the spectral shapes of 
the GVKM and Livermore fluxes presented in Section 2. 

6.2 Measurement of the neutron star mass 

A large set of neutron star masses have been previously measured and are listed in Ref. [40]. 
These measurements provide relative uncertainties at the 0.1% level on the neutron star mass 
for binary systems. However, an isolated system such as a single supernova without companion 
would not provide as much precision. Hence, our interest to retrieve the neutron star mass 
from the number of detected neutrino events. We now fix the supernova distance, assumed to 
be known by another conventional technique and try to retrieve the neutron star mass from 
the number of detected events and the expression of the neutron star gravitational energy. In 
a practical case, the measurement of the distance could be achieved through the Expanding 
Photosphere Method [42], This technique led to a 10 % uncertainty of the SN1987A distance. 
Using this as a reference, we fix the uncertainty on the reconstructed distance to 10%. We 
also fix the neutron star radius and associated uncertainty, as in Section 6.1. We can thus 
assess how well one could measure the neutron star mass by fitting the energy distribution 
of incoming neutrino events as in Section 6.1. We generated neutrino supernova spectra with 
statistics corresponding to different masses between 1 and 2 solar masses. For a canonical 
supernova located at 10 kpc, the mass could be measured with an uncertainty of 16% and 
13% for a 1 kt and a 20 kt scintillator detector respectively. 

2 In Sections 6.1 and 6.2 an approximated relation of the gravitational binding energy of a sphere of uniform 
density is employed. We use it here for illustrative purposes, since the determination of the supernova distance 
and of the neutron star mass are not the main purposes of the present work. 
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6.3 An example: Disentangling different supernova flux shapes 

The far location of SN1987A in the Large Magellanic Cloud has provided information on the 
total gravitational energy radiated as neutrinos and on the neutrino temperatures, consistent 
with predictions (see Ref. [4] and references herein). However, the statistic was not large 
enough to precisely determine the supernova neutrino energy spectrum. In what follows we 
study the possibility to discriminate between two flux shapes, GVKM and Livermore, with 
the current and future generation of LLSDs. We statistically compare spectra generated from 
a GVKM distribution with both GVKM and Livermore fluxes. To simplify the problem as 
well as remain conservative, we normalize both fluxes to the same fluence, only their spectral 
shapes differ. We then compare the y 2 distributions of 5,000 simulated experiments fitted 
with the GVKM and Livermore fluxes. The two y 2 distributions are generated for supernova 
distances ranging from 0 to 10 kpc. The separation power between the two is expressed 
through the following figure of merit: F.o.M. = ., 35 \x^+aTj w ^h /H ,2 and the respective 
mean and standard deviation of the two distributions. This figure of merit has been carefully 
chosen for its efficiency in disentangling Gaussian-like probability distributions. It has been 
used in the Nucifer experiment to separate neutron and gamma backgrounds via the Pulse 
Shape Discrimination (PSD) technique. More details can be found in Ref. [43]. Figure 14 
shows the behavior of this separation power parameter as a function of the increasing distance 
for a 1 kt and a 20 kt detector. The reconstruction of a GVKM spectrum at 95% C.L. can 
be achieved up to a distance of 2.2 kpc for a 1 kt detector such as KamLAND and 9.5 kpc 
for a 20 kt detector such as JUNO. 




Figure 14. Figure of merit for disentangling GVKM and Livermore energy distributions as a function 
of the supernova distance for a 1 kt (left) and a 20 kt (right) detector. The dotted line represents the 
upper limit over which a GVKM spectrum is reconstructed with at least 95% C.L. 


7 Discussion 

In this paper, we discussed the interest of combining worldwide large liquid scintillator detec¬ 
tors to efficiently reconstruct a supernova location. Such an information would be precious 
to the SNEWS network in order to provide not only an alert that the optical observation of 
a supernova is imminent but also that it is likely to occur in a restricted region of the sky. 
While for a typical supernova located at 10 kpc, a 45°(68% C.L.) reconstruction cone would 
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be currently too large to be entirely enclosed within the aperture of an optical telescope 3 , it 
would focus astronomers attention to a particular held of view, enhancing the probability for 
a comprehensive optical observation. The impending Zwicky Transient Facility (ZTF) [44] 
will be able to completely scan this region of the sky within two hours. The presence of 
several telescopes based on the design described in Ref. [45] would provide a high probability 
of successful observation as well. The implementation of future detectors would soon reduce 
the region of interest to a 12°(68% C.L.) radius cone. In the case of a nearby supernova 2 kpc 
away, a single detector such as KamLAND could point its direction within a 10°(68% C.L.) 
cone and the combinations of current and/or future detectors could reduce this uncertainty 
below 5°(68% C.L.). 

Besides providing a relevant direction reconstruction, a wealth of information can be 
obtained through IBD events. In particular, in Section 6 we have discussed the extraction 
of information on the supernova distance and the mass of the neutron star. We also dis¬ 
cussed how accurately two flux models could be discriminated, as a function of the supernova 
distance. 

Though it has hardly been discussed in this paper, the observation of a large quantity of 
neutrinos in several large detectors will help improve the direction reconstruction via triangu¬ 
lation technique. The combined data from several time-calibrated detectors scattered over the 
globe would be contained within the 40 ms time window corresponding to the Earth’s cross¬ 
ing. The detection of a full neutrino spectrum at each moment of the supernova luminosity 
rise and decay would decrease any statistical uncertainty related to the luminosity rise time 
thus making triangulation a viable technique as explained in Ref. [32], Since only IceCube, 
Super-Kamiokande and NOzzA’s far detector are currently operating with active masses rel¬ 
evant for triangulation, it might not be a suitable method as we speak. However, with the 
impending construction of large detectors such as JUNO, LENA or Hyper-Kamiokande, it 
could be considered as reliable as the other methods. Another pointing method relying on the 
observation of spectral distortions caused by matter effects in Earth is discussed in Ref. [46]. 

Recent technological developments have been successfully carried out in the held of 
detection materials. For example, detectors using liquid scintillators mixed and dissolved into 
water are currently being designed. A detailed study of such technology and its potential 
physics applications can be found in Ref. [47]. Besides combining the technical advantages 
of both materials, a water-based liquid scintillator detector also provides two directional 
information. Indeed, enhanced by the possible use of fast photodetectors such as LAPPD’s 
(Large Array Picosecond PhotoDetectors), the collection of Cerenkov and scintillation light 
produced by the both IBD and electron scattering processes would exploit the full directional 
information carried by neutrinos. 

We would like to draw attention to technical issues that might arise along with the high 
event rates in large detectors. Underground detectors are designed to operate in an ultra low 
background environment with typical trigger rates of a few hundred Hz. The detection of tens 
of thousands of neutrinos even distributed in time over a 10 s window with a strong peak at 
1 s may lead to a non-negligible dead-time or, worst case scenario, to a DAQ crash. On the 
other hand, the NOz'A far detector, since it operates at the surface, is subjected to a 300 kHz 
background rate from cosmic rays at 10 MeV. While thanks to its fine-grained segmentation, 
it could separate IBD events from cosmic muons tracks, an efficient supernova tagging would 
place significant demands on the DAQ in terms of dead-time and data readout and analysis 

3 Note that the Moon’s angular diameter over the sky is 0.5° hence enclosed in a 0.25° radius cone. 
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rates. In order for any of these detectors to efficiently collect a large number of neutrinos in 
a short period of time and increase its pointing ability, considerable efforts will be required 
to improve the current state of data acquisition. 

8 Conclusion 

In conclusion simultaneous detection of several hundreds or thousands of neutrinos via several 
processes will bring an enormous amount of information about stellar evolution and neutrino 
flavor conversion occurring in extreme media. The directional information extracted from 
these same neutrinos is a precious input for astronomers and could lead to a prompt and 
therefore comprehensive observation of a galactic supernova in several domains of the electro¬ 
magnetic spectrum. To link these two observations and recover as much data as possible from 
the next and long-awaited supernova in our galaxy, the SNEWS network performs a constant 
watch that could be improved by incorporating all large neutrino detectors, all technologies 
combined, present or future. 
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Appendices 

A Toy Monte-Carlo 

In order to generate an important amount of supernova neutrino events, we oriented our choice 
toward a fast toy Monte-Carlo simulation rather than a full GEANT4 [48]. It implies a good 
understanding of IBD kinematics and detection to only keep its most important tendencies, 
without losing too much information. The code is written in MATLAB and uses GEANT4 
outputs. 

A.l Event generator 

Our toy Monte Carlo generates a set of neutrinos picked in a supernova spectrum flux. The 
(neutron, positron) pairs momenta and directions are created according to kinematics exposed 
in Section 4.1. Their respective detection vertices are generated using data from Figures 6 
and 8. In the end, our toy Monte Carlo supplies a visible energy as well as prompt and 
delayed vertices. 

Neutrino generation 

Neutrinos energies E u are picked at random in a supernova spectrum among the ones presented 
in Section 2. Depending on the analysis, the number of generated events can either be user- 
defined or computed with respect to a defined supernova distance and detector size. 
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(e + ,n) pair generation 

For each event, a positron angle is picked using probability density function defined in Sec¬ 
tion 4.1. From this angle and the initial E u , the positron momentum p e = \JE% — is 
defined using Eq. 3.2. Considering protons to be at rest in the laboratory frame, the neutron 
momentum is computed through the p n = pp e — p e relation. Until now, the generation of the 
(e + ,n) pair was purely dictated by analytical IBD kinematics. However in order to obtain 
realistic prompt and delayed events, it is mandatory to simulate the behavior of the positron 
and neutron in a liquid scintillator. For simplicity reasons and since it is the current detector 
expected to detect the most neutrino events from a nearby supernova with the best sensitiv¬ 
ity, we carried out simulations using the KamLAND liquid scintillator composed of a mixture 
of 80% of C 12 H 20 and 20% of C 9 H 12 with a density of 0.78 g/crn 3 . After cross-checking the 
results using another liquid scintillator, we noticed no significant discrepancy between the 
two liquids and therefore used the KamLAND liquid scintillator for all simulations. 

Prompt event 

Given their ionizing and penetrating nature, the path of positrons in the toy Monte Carlo 
is considered a straight line whose direction is given by the initial positron angle picked 
during the IBD generation. By simulating positrons with energies between 0 and 100 MeV 
in GEANT4, we obtained the mean positron range distribution with respect to the positron 
energy. In the toy Monte Carlo, the prompt event vertex is located on the positron direction 
and its distance from the interaction point is picked at random in the range distribution 
corresponding to its energy. Note that after it deposits its energy along its pathway, the 
positron annihilates with an electron thus creating 2 gammas of total energy 2 m e , with m e 
the electron mass. Hence the energy deposition of the positron is not point-like but has a 
shape that depends on the initial positron energy. However, liquid scintillator detectors cannot 
resolve this shape and the reconstructed vertex usually corresponds to the barycenter of the 
charges. To efficiently take this effect into account in the toy Monte Carlo, we consider the 
path length of positrons from their creation to their annihilation and then apply a correction 
factor on the track accounting for point-like reconstruction of the positron track. Using 
GEANT4, the value of this correction factor was found to be equal to 0.72. The reconstructed 
prompt vertex is then smeared accordingly to the user-defined vertex resolution. 

Delayed event 

As stated in Section 4.1, the neutron loses its energy via elastically scattering off the nuclei 
composing the liquid, mainly light H nuclei. Given the small mass difference between a 
neutron and a H nucleus, neutrons lose half their kinetic energy, on average, at each scattering. 
Neutrons become rather quickly thermalized (T n ~ 0.025 eV) and each scattering saves some 
directional information of the initial momentum, the first one being the most significant. 
After being thermalized, the neutron diffuses in the liquid before it gets captured on a H or 
Gd atom. This capture leads to the emission of gamma rays whose path distribution needs to 
be taken into account as well. In our GEANT4 simulations, we first obtain the path length 
until thermalization for neutrons of different energies as displayed in Figure 8 (left). Then, by 
shooting thermal neutrons from this point, we obtain the diffusion path length before capture 
in an unloaded or a Gd-loaded liquid as seen in Figure 8 (right). Finally, in the case of a 
capture on H, we shot 2.2 MeV gamma rays to obtain their path length distribution. In the toy 
MC, the thermalization vertex is generated with respect to the mean path length depending 
on the neutron initial energy. The delayed vertex is the combination of the thermalization 
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vertex and two isotropic vectors whose respective lengths are picked from the thermal neutron 
diffusion distribution and the gamma ray path length distribution. Like the prompt event, 
the reconstructed delayed vertex is smeared by the vertex resolution. 



Figure 15. Mean projections of prompt and delayed events along the v e direction obtained with 
a GEANT4 simulation and our tuned toy MC. The projections are expressed as a function of the 
neutrino energy. 


A.2 Toy MC reliability 

To assess our toy MC IBD detection outline reliability, comparisons are made with GEANT4, 
for uniformly distributed neutrinos between 0 and 100 MeV, relevant energy for supernova 
neutrinos. 

The most relevant observables that can be used to compare the Delayed-Prompt vector 
in both GEANT4 and our toy MC are the prompt and delayed projections of the incoming 
neutrino axis. A 20 cm vertex resolution has been applied to all simulations. Figure 15 shows 
the good agreement between our two simulation methods. 
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